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ABSTRACT 

Low-mass planets are known to undergo Type I migration and this process must have played 
a key role during the evolution of planetary systems. Analytical formulae for the disc torque 
have been derived assuming that the planet evolves on a fixed circular orbit. However, recent 
work has shown that in isothermal discs, a migrating protoplanet may also experience dynam¬ 
ical corotation torques that scale with the planet drift rate. The aim of this study is to examine 
whether dynamical corotation torques can also affect the migration of low-mass planets in 
non-isothermal discs. We performed 2D radiative hydrodynamical simulations to examine the 
orbital evolution outcome of migrating protoplanets as a function of disc mass. We find that 
a protoplanet can enter a fast migration regime when it migrates in the direction set by the 
entropy-related horseshoe drag and when the Toomre stability parameter is less than a thresh¬ 
old value below which the horseshoe region contracts into a tadpole-like region. In that case, 
an underdense trapped region appears near the planet, with an entropy excess compared to 
the ambient disc. If the viscosity and thermal diffusivity are small enough so that the entropy 
excess is conserved during migration, the planet then experiences strong corotation torques 
arising from the material flowing across the planet orbit. During fast migration, we observe 
that a protoplanet can pass through the zero-torque line predicted by static torques. We also 
hnd that fast migration may help in disrupting the mean-motion resonances that are formed 
by convergent migration of embryos. 

Key words: accretion, accretion discs - planet-disc interactions- planets and satellites; for¬ 
mation - hydrodynamics - methods: numerical 


1 INTRODUCTION 


A striking feature of the population of exoplanets discovered so far 
is their great diversity. Compact and non-resonant systems of super- 
Earths and mini-Neptunes have heen discovered orbiting within a 
few tenths of AU from their stars (e.g. Lissauer et al. 2011; Lovis et 
al. 2011). Giant planets with periods from days to several thousands 
of days have also heen found around 14% of Sun-like stars, with 
a frequency that is strongly correlated to the metallicity of the host 
star (Cumming et al. 2008; Mayor et al. 2011). A main aim of planet 
formation and evolution scenarii is to examine whether this broad 
diversity can be explained, and under which conditions. These sce¬ 
narii can be constrained using for example planet population syn¬ 
theses which combine analytical models of disc evolution, planet 
formation and orbital evolution due to disc-planet interactions (Ida 
& Lin 2008; Mordasini et al. 2009; Mordasini et al. 2012); and gen¬ 
erate a large variety of system architectures that can be statistically 
compared to observations. An alternative method is to make use of 
N-body simulations of planetary systems formation which do not 
cover a parameter space as large as in planet population syntheses, 
but allow for a self-consistent treatment of planet-planet interac¬ 


tions (Hellary & Nelson 2012; Coleman & Nelson 2014; Cossou et 
al. 2014). 

In both of these approaches, disc-planet interactions causing 
planet migration are modeled using analytical formulae for the disc 
torque experienced by the planet. This gravitational torque consists 
of two components. The differential Lindblad torque results from 
the angular momentum exchange between the planet and the spi¬ 
ral density waves it generates inside the disc. For sufficiently low- 
mass planets, it scales linearly with disc mass and planet mass and 
as the inverse square of the disc aspect ratio (Tanaka et al. 2002). 
Although its sign depends on the density and temperature gradi¬ 
ents inside the disc, the differential Lindblad torque is generally 
negative for typical disc models and is therefore responsible for 
inward migration. The corotation torque is due to the torque ex¬ 
erted by the material located in the coorbital region of the planet. 
It is composed of a barotropic part which scales with the vorten- 
sity (i.e. the ratio between the vertical component of the vortic- 
ity and the disc surface density) gradient (Goldreich & Tremaine 
1979) plus an entropy-related part which scales with the entropy 
gradient (Baruteau & Masset 2008; Paardekooper & Papaloizou 
2008). A negative vortensity (resp. entropy) gradient gives rise to 
a positive vortensity (resp. entropy) related corotation torque. It 
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has been shown that for midly positive surface density gradients or 
negative entropy gradients, a positive corotation torque can eventu¬ 
ally counteract the effect of a negative differential Lindblad torque, 
which may stall or even reverse migration (Masset et al. 2006; 
Paardekooper & Papaloizou 2009). In isothermal discs, the corota¬ 
tion torque is a non-linear process generally referred as the horse- 
hoe drag and whose amplitude is controlled by advection and dif¬ 
fusion of vortensity inside the horsehoe region. In non-isothermal 
discs, the corotation torque is also powered by singular produc¬ 
tion of vortensity due to an entropy discontinuity on downstream 
separatrices (Masset & Casoli 2009; Paardekooper et al. 2010). In 
the absence of any diffusion processes inside the disc, vortensity 
and entropy gradients across the horseshoe region tend to flatten 
through phase mixing, which causes the two components of the 
horseshoe drag to saturate. Consequently, desaturating the horse¬ 
shoe drag requires that some amount of viscous and thermal dif¬ 
fusions are operating inside the horseshoe region. In that case, 
the amplitude of the horseshoe drag depends on the ratio between 
the diffusion timescales and the horseshoe libration timescale and 
its optimal value, also referred as the fully unsaturated horseshoe 
drag, is obtained when the diffusion timescales are approximately 
equal to half the horseshoe libration time (e.g. Baruteau & Masset 
2013). In the limit where the diffusion fimescales become shorter 
than the U-tum timescale, the corotation torque decreases and ap¬ 
proaches the value predicted by linear theory. Therefore, the coro¬ 
tation torque can be considered as a linear combination of the fully 
unsaturated horseshoe drag and the linear corotation torque with 
coefficients depending on the ratio between the diffusion timescales 
and the horseshoe libration timescale. Corotation torque formulae 
as a function of viscosity and thermal diffusivity were recently pro¬ 
posed by Paardekooper et al. (2011) and Masset & Casoli (2010). 

Such analytical formulae have been derived assuming that the 
planet evolves on a fixed circular orbit. In that case, the planet feels 
a static torque that depends only on the local temperature and sur¬ 
face density conditions in the disc. If the planet is allowed to mi¬ 
grate, however, the disc torque can also exhibit, under certain con¬ 
ditions, a dependence on the planet drift rate. For instance, this can 
occur if the planet is massive enough to partly deplete its coorbital 
region. In that case, a coorbital mass deficit is created and the coro¬ 
tation torque is essentially due to the gas material that flows across 
fhe planet’s horseshoe region as the latter migrates. This contribu¬ 
tion to the corotation torque not only scales with the planet drift 
rate, but also has the same sign as the drift rate. This yields a pos¬ 
itive feedback on migration, and eventually gives rise to runaway 
effects if the disc is massive enough (Masset & Papaloizou 2003). 

In a recent study, Paardekooper (2014) has shown that low- 
mass planets embedded in isothermal discs a few times more mas¬ 
sive than the Minimum Mass Solar Nebula (hereafter MMSN; 
Hayashi 1981) may also experience runaway migration in some 
cases. For low-mass planets that do not significantly perturb the 
underlying disc structure, this process does not result from a coor¬ 
bital mass deficit but rather arises due to the action of dynamical 
corotation torques that scale with the planet drift rate and depend 
on the vortensity gradient inside the disc (Paardekooper 2014). In 
the case where the migration proceeds in the direction set by the 
corotation torque, Paardekooper (2014) has shown that these dy¬ 
namical torques have a positive feedback on migration, and can 
possibly give rise to high drift rates. Interestingly, using an isother¬ 
mal disc model with a positive surface density gradient to mimic 
what may happen in non-isothermal disc, Paardekooper (2014) has 
shown that an outward migrating protoplanet subject to strong dy¬ 
namical torques can pass through the location of the zero-torque 


radius where the Lindblad torque and the (static) corotation torque 
cancel each other. 

In this paper, we focus on the case of non-isothermal disc 
models, in which outward migration can proceed due to the opera¬ 
tion of the entropy-related horseshoe drag. The aim of this study is 
to examine the plausibility that low-mass planets embedded in non- 
isothermal discs experience strong dynamical corotation torques as 
they migrate. Using hydrodynamical simulations, we find fhat in 
relafively massive discs with an initial large-scale entropy gradient, 
the drift rates reached by a migrating planet can indeed be signifi¬ 
cantly higher than those expected by assuming a static torque. We 
show that this occurs when i) the planet migrates in the direction set 
by the entropy-related horseshoe drag, and ii) the Toomre stability 
parameter is less than a critical value that depends on the entropy 
gradient. Below this value, the horseshoe region does not extend to 
the full In in azimuth and an underdense trapped region appears 
which causes the planet to undergo large drift rates. In this fast mi¬ 
gration regime, we also find thaf an oufward migrafing protoplanef 
can migrate well beyond the location of the zero-torque radius, sim¬ 
ilarly to the findings of Paardekooper (2014) in isothermal discs. 

This paper is organized as follows. In Sect. 1, we present the 
numerical setup. In Sect. 2, we present the results of the simulations 
for non-isothermal disc models. The case of radiative disc models 
that include the effects of viscous plus stellar heating, and radiative 
cooling are discussed in Sect. 3. Finally, we discuss our results and 
draw our conclusions in Sect. 4. 


2 THE HYDRODYNAMICAL MODEL 

Simulations were performed using the GENESIS (De Val-Borro 
et al. 2006) numerical code which solves the equations governing 
the disc evolution on a polar grid {R, tp) using an advection scheme 
based on the monotonic transport algorithm (Van Leer 1977). It 
uses the FARGO algorithm (Masset 2000) to avoid time step limi¬ 
tation due to the Keplerian velocity at the inner edge of the disc. 

We adopt computational units such that the mass of the central 
star is M, = IMq, the gravitational constant is G = 1, and the radius 
R = 1 in the computational domain corresponds to 5.2 AU. When 
discussing the results of the simulations, time will be measured in 
units of the orbital period at R = 1. 

We do not consider the disc self-gravity in this work. The typ¬ 
ical value for the Tooomre stability parameter Q = kc JUG'S., where 
K is the epicyclic frequency, Cj the sound speed and S the surface 
density is G 3-4, suggesting thereby that the effect of self¬ 
gravity may be important. The gravitational potential for the disc 
therefore restricts to the gravitational potential from the star and 
the planet, plus the indirect terms that account for the fact that the 
frame centered on the star is not inertial. The gravitational poten¬ 
tial of the planet is smoothed over a distance b = QAHp, where 
Hp is the disc scale height at the location of the planet. Here, the 
planet orbit can evolve due to the gravitational forces exerted by the 
star and the disc, and eventually due to the gravitational interaction 
of the planet with other bodies (see Sect. 14.2t . When calculating 
the disc force experienced by the planet, we follow Paardekooper 
(2014) and exclude the axisymmetric component of the force to 
ensure consistency as the planet migrates. 

In the following, two different forms of the energy equation 
will be employed. We will first consider non-isothermal disc mod¬ 
els for which the energy equation corresponds to an advection- 
diffusion equation for the gas entropy with a constant thermal dif¬ 
fusion coefficient (Paardekooper et al. 2011). Then we will focus 
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on more realistic, radiative disc models whose thermal structure 
results from the balance between viscous plus stellar heating, and 
radiative cooling. 

2.1 Non-isothermal discs 

2.1.1 Numerical setup 

For non-isothermal disc models, we solve the following energy 
equation: 

de , 

— + V • (ev) =-(y-l)eV • v-H^eV-log5 (1) 

at 

where e is the thermal energy density, v the velocity, y the adiabatic 
index which is set to 7 = 1.4, and where ^ is a thermal diffusion 
coefficient. In the previous equation, S = p/l/ is the gas entropy, 
with p = (y - l)e the gas pressure. Strictly speaking, the diffusion 
term in the previous equation should involve the laplacian of tem¬ 
perature rather than entropy. This form for the energy equation is 
used so that it corresponds to an advection-diffusion equation for 
the gas entropy, and is valid provided that the variations in temper¬ 
ature occur over a scale smaller than those corresponding to the gas 
pressure (Masset & Casoli 2010). Here, this condition if fulfilled 
since we are interested in low-mass planets for which the width 
of the horseshoe region is typically a fraction of the disc pressure 
scale height. The thermal diffusion coefficient is chosen assum¬ 
ing a Prandtl number Pr = v/x of unity, where v is the kinematic 
viscosity which is assumed to be constant. In the case where dif¬ 
fusion of entropy is due to turbulent transport, this is consistent 
with the results of Pierens et al. (2012) who found Pr ~ 1.2 in 
non-isothermal disc models with turbulence driven by stochastic 
forcing. In a real protoplanetary disc, however, where radiative dif¬ 
fusion can possibly dominate over turbulent transport for diffusing 
entropy, we expect the Prandtl number to be of the order of unity 
only in the optically thick inner regions (Pierens et al. 2012), while 
a value < 1 is expected in the outer, optically thin regions. 

Typically Nr = 1600 radial grid cells uniformly distributed 
between P,,, = 0.4 and Ra„, = 5 are used, and = 2130 azimuthal 
grid cells. 

We employ closed boundary conditions at both the inner and 
outer edges of the computational domain, and make use of wave¬ 
killing zones for R < 0.5 and R > 4.5 to avoid wave reflections at 
the disc edges. 


2.1.2 Initial conditions for the non-isothermal runs 

The initial disc surface density is S = Eo(P/ao)“'^, where Eq is 
the surface density at the initial position of the planet R = ao and 
where the power law index is such that cr = 3/2 by default. For 
such a disc surface density profile, the initial vortensity-related part 
of the corotation torque cancels out so that the corotation torque 
consists only of its entropy-related part. Following Paardekooper 
(2014), we parametrise the disc mass through the parameter 


for which we adopt values ranging from qj = 0.005 to = 0.03. 

The initial aspect ratio is h = hoiR/ao)! where hg = 0.05 and 
/ is the disc flaring index for which we consider values of / = 0.7 
and / = -0.8. This corresponds to initial temperature profiles that 
vary as R^^ with p = 1 - 2/ = -0.4 and f = 2.6 respectively, and to 
initial entropy profiles that vary as Rr^ with f = p — (y — Vja = 
and ^ = 2 respectively. 


Our reference value for the planet-to-star mass ratio is g = 
10“^ (equivalent to a 3.3 Earth mass planet in our units). For ho = 
0.05, the half-width of the horseshoe region is ~ l.lflo - fqjho ~ 
0.016 ao (Masset et al. 2006), which means that the horseshoe re¬ 
gion is typically resolved by about 12 grid cells in the radial direc¬ 
tion. 


2.2 The radiative disc model 

2.2.1 Numerical setup 


In addition to these non-isothermal disc models, we have also con¬ 
sidered more realistic discs in which the effects of viscous and stel¬ 
lar heating are balanced by radiative cooling. In that case, the fol¬ 
lowing energy equation is adopted: 

^ + V . (cv) = -(7 - l)eV • V -H Qf -Q--2HW F (2) 

at 

where H is the disc scale height. In the previous equation, is the 
viscous heating term, = IcrsTf^ is the local radiative cooling 
from the disc surfaces, where ctr is the Stephan-Boltzmann con¬ 
stant and Teff the effective temperature which is given by (Menou 
& Goodman 2004): 


pA _p4 V 3 1 

rfr-^ wtth r.,y--rv-.- (3) 


Here, T is the midplane temperature and t = kI,I 2 is the vertical 
optical depth, where k is the Rosseland mean opacity which is taken 
from Bell & Lin (1994). r,>r is the irradiation temperature which is 
computed from the irradiation flux (Menou & Goodman 2004): 


<^bTI=A 


4nR^ R 


I dlogH \ 
[dlogR j 


(4) 


where e = 1/2 is fhe disc albedo and L* fhe stellar luminos¬ 
ity which is computed assuming a stellar radius R* = 1.5 Rq 
and a stellar temperature T* = 4370 K. The factor A is set to 
A = 1 if d(HIR)ldR > 0, whereas we set A = 0 in the regions 
where d(HIR)ldR < 0 in order to take into account possible self¬ 
shadowing effects (Gunther & Kley 2004). In Eq.|2] F is the radia¬ 
tive flux which is treated in the flux-limited diffusion approach and 
which reads (e.g. Kley & Crida 2008): 


F 


l6crAF 

-vr 

pK 


(5) 


where p = E/2// is the mid plane density and where T is a flux- 
limiter (e.g. Kley 1989). 


2.2.2 Initial conditions for the radiative runs 

For the radiative simulations, we will focus on the case of equilib¬ 
rium discs with surface density E = Eo(////?o) and constant vis¬ 
cosity, where Eo is the surface density atRo = 1. Two disc models 
are considered, with parameters summarized in Table 1. For Model 
1, Eq = 1.6 X 10“^ (equivalent to Eq = 520 g.cm^^ at 5.2 AU) and 
V = 10“^, whereas for Model 2, Eq = 3.2 x 10“^ (equivalent to 
Eo = 1024 g.cm-2 at 5.2 AU) and v = 10^®. 

The resolution used in each simulation is {Nr,NP = 
(1024,1536) with a radial extent corresponding to R e [//,„,//„„;], 
where /?,„ and R„,„ are given in Table 1. 

We adopt planet masses such that q = 3x 10“^, 6 x 10“^, 10“"* 
for Model 1 andq = 10“^, 1.5x 10“^,2x 10“^ for Model 2. For each 
model, the planet’s horseshoe region is resolved at minimum by 
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Figure 1. Left panel: Time evolution of the semi-major axis (top) and drift rate as a function of radius (bottom) of a protoplanet with q = 10“^ embedded in a 
non-isothermal disc with ^ = —1 and v = x = 10“^, for various disc masses. The dashed line shows the drift rate that is expected from analytical formulae for 
the disc torque, whereas the dot-dashed line corresponds to the drift rate expected assuming a fully unsaturated entropy-related horseshoe drag. Right panel: 
same but for an initial entropy gradient corresponding to f = 2. 


Table 1. Parameters of the radiative simulations, where Z = I.o{R/Ro) 
and where the viscosity is constant 


Model 

So 

V 

Rin 

^OUT 

Model 1 

1.6 X 10^3 

10-3 

0.4 

6 

Model 2 

3.2 X 10-3 

10-® 

0.4 

3.5 


~ 11 grid cells along the radial direction, which leads to a relative 
error on the estimation of the corotation torque of ~ 10% (Masset 
2002 ). 

3 NON-ISOTHERMAL DISC MODELS 
3.1 Results 

For a surface density profile with cr = 1.5, the total disc torque F 
exerted on a protoplanet held on a fixed circular orbit, which in the 
following we will refer to as the static torque, can be considered as 
the sum of the differential Lindblad torque, and the entropy-related 
corotation torque. In the limit where the corotation torque is fully 
unsaturated, F is given by (Paardekooper et al. 2010): 

yF/Fo =-2.5-1.7/3 -fO.lcr-H 7.9-, (6) 

r 

with Ffl = {qlhpy-'Lpa^Q.^p. Here a is the planet semi-major axis and 
the subscript ”p” indicates that quantities are evaluated at the or¬ 
bital position of the planet. For ^ = -1, this gives yF/Fo-7.3, 


whereas yF/Fo ~ 4.5 for ^ = 2. For these two disc models, the to¬ 
tal and entropy-related corotation torques have therefore the same 
sign, and the aim of this section is to examine whether in that 
case, a migrating protoplanet can reach migration rates that are 
significantly higher than those obtained by evaluating the static 
torque only. In isothermal discs a few times more massive than the 
MMSN, it has been shown that if the planet migrates in the direc¬ 
tion set by the corotation torque, it can be subject to strong dynam¬ 
ical corotation torques and runaway effects can even be triggered 
under certain conditions (Paardekooper 2014). 

The upper left panel of Fig. [T] shows, for the disc model with 
^ = -I and for various disc masses, the semi-major axis evolu¬ 
tion of a ^ = 10“^ planet initially located at qq = 3. In agree¬ 
ment with the above estimation of the total torque, the planet mi¬ 
grates inward in that case, and in the direction set by the static 
corotation torque since the entropy gradient is positive. Here, both 
the viscosity and thermal diffusivity are set to = 10 “®, consis¬ 
tently with the assumption of a Prandtl number of unity. This value 
for the thermal diffusivity is chosen such that the entropy-related 
horseshoe drag is initially close to its unsaturated value. This is 
reached when the thermal diffusion timescale across the horse¬ 
shoe region approximately equal to half the libration 

timescale Tm, = Snal(3Q.pXs). However, it can be easily shown that 
TdlTiib oc which evaluates to Tj/T/,i, oc for / = 0.7. 

Therefore, we expect the level of saturation of the corotation torque 
to increase as the planet migrates inward. In the lower left panel of 
Fig.m are compared, for the runs with = 0.005,0.02, the mi¬ 
gration rates at different locations with those predicted from static 



















Fast migration of low-mass planets in radiative discs 5 


torques. Also plotted on this figure are the drift rates obtained as¬ 
suming a fully unsaturated horseshoe drag, and which should there¬ 
fore correspond, in principle, to the highest drift rates that could be 
reached by a migrating protoplanet. Clearly, the migration rates are 
significantly higher in the case where the planet is allowed to mi¬ 
grate, with values that appear to be even larger than the estimate 
given by assuming a fully unsaturated corotation torque. 

The right panel of Fig.[T]presents the results of the simulations 
for the disc model with ^ = 2 and for which outward migration is 
expected. We note that the run with = 0.03 serves for illustrative 
purpose only since < 1 for > 2 in that case. Here oc a'-’, 

leading to an increase in the saturation of the corotation torque as 
the planet migrates outward. The main implication of this is that 
the total torque acting on the planet can eventually cancel at the 
location where the differential Lindblad torque and the saturated 
corotation torque have the same amplitude. An estimation for the 
location of this zero-torque radius as can be obtained by equating 
the diffusion and libration timescales across the horseshoe region. 
It is straightforward to show that this gives: 


flo 



(7) 


where Q.o is the angular velocity at the initial location of the planet. 
For;^' = 10“*, q = 10“*, and / = -0.8 we find aj/flo ~ 1.6. Al¬ 
though not shown here, we find that such an estimation agrees fairly 
well with the value for the zero-torque radius deduced from runs 
with the planet held on a fixed orbit. The final evolution outcome 
for runs with qj < 0.01 remains uncertain, but we can reasonably 
expect that the migration will ultimately be halted at a ~ as for 
these cases. For q^ > 0.02, however, the migration rate continu¬ 
ously increases as the evolution proceeds, leading to the planet mi¬ 
grating outward very rapidly. This makes the planet pass through 
the location of the zero-torque radius predicted by static torques 
and rapidly reach the outer edge of the computational domain. In 
that case, the migration rate is again higher than what expected 
from the action of the Lindblad plus the fully unsaturated entropy- 
related corotation torque, as illustrated by the lower right panel of 
Fig.[T]which shows the drift rate as a function of radius for the sim¬ 
ulation with qj = 0.02. This result seems to indicate that in non- 
isothermal discs and provided the disc mass is high enough, a pro¬ 
toplanet can feel additional dynamical corotation torques as it mi¬ 
grates, similarly to what occurs in isothermal discs (Paardekooper 
2014). 

In Fig. 12 is presented the orbital evolution of a ^ = 3.10“* 
planet for the same disc model, except that both the viscosity and 
thermal diffusivity have been increased to v = ^ = 10“* to make 
sure that the underlying surface density profile is not altered by 
the presence of the planet. Eq. |3 predicts that outward migration 
should come to a halt at as ~ 2.35 in that case. Although it seems 
reasonable to claim that migration will stall close to this location in 
the simulation with q^ = 0.005, we see that the planet undergoes 
fast outward migration and crosses the location of the zero-torque 
radius in runs with > 0.01. Here ~ 0.8, so that it is plausible 
that here, the occurence of fast outward migration for lower disc 
masses might be related to the onset of non-linear effects, resulting 
in a boost of the corotation torque (Masset et al. 2006). 


3.2 Coorbital mass deficit 

In order to get some insight about the origin of the observed large 
drift rates, we display in Fig.[2the perturbed surface density near 



Figure 2. Time evolution of the semi-major axis of a planet with q = 10“* 
embedded in a non-isothermal disc with f = 2 and v = x = 10“*, for various 
disc masses. 


the planet at r = 200 orbits, for the disc model with ^ = -1 (left 
panel) and f = 2 (right panel). For^ = -1 (resp. ^ = 2) the initial 
entropy gradient is positive (resp. negative) and dynamics in the 
horseshoe region yields a positive (resp. negative) surface density 
perturbation at the outward downstream separatrix whereas a neg¬ 
ative (resp. positive) surface density perturbation is created at the 
inward downstream separatrix. Overplotted as black lines are a few 
streamlines which correspond to material that flows across the or¬ 
bit by executing a single U-turn, whereas the overplotted white line 
shows the librating material bound to the planet. We see that the 
radial drift of the planet causes the horseshoe streamlines to be sig¬ 
nificantly distorted and to contract into a tadpole-like region. This 
arises because the planet drift rate a is higher than the critical drift 
rate df above which the planet migrates over a distance greater than 
the horseshoe width over one horseshoe libration time. 

As the planet migrates, the material that flows across the or¬ 
bit has always a positive feedback on migration, whereas librating 
material that moves with the planet has always a negative feedback 
on migration. However, in the case where the planet migrates in the 
direction set by the corotation torque, it is clear from Fig. [2 that 
the tadpole-like region tends to be underdense compared to the rest 
of the disc, resulting in a net positive feedback on migration, and 
possibly leading to high drift rates. We note that this effect cannot 
occur if the planet migrates in the opposite direction to that set by 
the corotation torque, since in that case the librating region would 
be over dense relative to the ambient disc, resulting in a significant 
negative feedback on migration. 

Interestingly, the presence of such an underdense librating re¬ 
gion at the leading side of the planet has also been observed in 
simulations of giant planets which undergo inward runaway mi¬ 
gration (Artymowicz 2004; D’angelo & Lubow 2008). Here, it is 
worth noting that the low-density region does not result from the 
gap opening process but is rather due to radiative effects. Never- 
thess, using the terminology employed in the framework of Type 
III migration. Fig. [2 shows that radiative effects can give rise to a 
coorbital mass deficit given by (Masset & Papaloizou 2003): 


6m = AnaXs (S^ - E*) (8) 

where E, is the surface density at the upstream separatrix and E;, the 
mean surface density in the librating region. In order to provide an 
estimation for the coorbital mass deficit in the librating region, we 
make in the following the assumption that the only contribution to 
the density perturbation in the librating region results from entropy 
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Figure 3. Left panel: Contours of the perturbed surface density at t = 200 orbits near the planet for the disc model with f = -1 and qj = 0.02. Overplotted as 
white lines are a few streamlines that correspond to the librating trapped material that moves with the planet and which has therefore a negative feedback on 
migration. The black lines are streamlines that correspond to the material that flows across the orbit as the planet migrates and which has a positive feedback 
on migration. In the case where the planet migrates in the direction set by the corotation torque, the trapped region is gas depleted, yielding a net positive 
feedback on migration. Right panel: same but for the disc model with f = 2 and qj = 0.02. 


advection. In reality, however, the density perturbation features an 
additional contribution linked to the production of vortensity at the 
outgoing separatrix (Paardekooper et al. 2011). Ignoring this con¬ 
tribution gives the following expression for the surface density E_ 
in the librating region (Baruteau & Masset 2008; Paardekooper et 
al. 2011): 

= E, |l - 2^-^-J—j for 0 < X < Xj (9) 

Here, we have also ignored the radial density gradient and assumed 
that the corotation torque was initially unsaturated. This gives a 
mean surface density in the librating region: 

I f Us 

E,. = E,(l - (10) 

ya 

The coorbital mass deficit then becomes: 

5m = An\^\x%ly ( 11 ) 

which shows that the coorbital mass deficit scales with the entropy 
gradient and the width of the horseshoe region. We note in passing 
that at the boundary corresponding to d ~ df, we expect the coro¬ 
tation torque exerted on the planet to be given by (e.g. Papaloizou 
et al. 2007): 

Tcr = ^dfilpOSm (12) 

Using the expression for the coorbital mass deficit given by Eg. IIII 


the previous equation then becomes: 

TCT^^^E.xtn" (13) 

4 y 

so that we recover the standard expression for the entropy-related 
horseshoe drag (Baruteau & Masset 2008). For d > df, however, 
the expression for the unsaturated corotation torque given by Eq. 
1121 is no longer valid and should be modified (Papaloizou et al. 
2007). In that case, we expect the corotation torque to be prevented 
from saturation since the horseshoe region does not longer extend 
to the full 2n in azimuth, so that phase-mixing cannot occur. In 
isothermal discs with a large-scale vortensity gradient, Ogilvie & 
Lubow (2006) found that drift rates with d > df can even lead 
to torques that are much larger than the unsaturated value in ab¬ 
sence of migration. They showed that this arises due to corotational 
torques caused by a vortensity asymmetry in the coorbital region 
between gas on the leading and trailing sides of the planet. In the 
limit of low viscosity, the trapped material tends to conserve its ini¬ 
tial vortensity as the planet migrates, leading to a possibly signifi¬ 
cant vortensity contrast between the librating region and the local 
disc. In that case, the impact on migration depends in fact on the 
sign of the vortensity deficit (Paardekooper 2014) which is defined 
as: 

(5(w/E) = 1 - (14) 

tuo/Eo 

where Wo/E„ is the local disc vortensity at the position of the planet 
and tuo/Eo the vortensity at the initial location of the planet. In the 
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case where the planet migrates in the direction set by the vortensity- 
related horseshoe drag, the vortensity deficit is positive and this 
yields a positive feedback on migration whereas if the planet mi¬ 
grates in the opposite direction to that set by the corotation torque, 
the vortensity deficit is negative and this gives rise to a negative 
feedback on migration (Paardekooper 2014). 

Here, the appearance of the streamlines in Fig. [3] combined 
with the very high drift rates (higher than those expected from 
unsaturated torques) that have been reported in the previous sec¬ 
tion, strongly suggests that a mecanism similar to that presented in 
Ogilvie & Lubow (2006) is at work here. However, since there is 
no initial vortensity gradient inside the disc, we suggest the strong 
corotation torques acting on the planet as being generated by an en¬ 
tropy rather than a vortensity deficit between the trapped material 
and the ambient disc. For non-isothermal discs, we therefore make 
use of an entropy deficit 6S that we define as; 

<55 = 1 - (15) 

where 5^ is the local disc entropy at the position of the planet and 
So the disc entropy at its initial location. In the case where the 
entropy deficit is positive, the planet migrates in the direction set 
by the entropy-related horseshoe drag and an underdense region 
with an entropy excess compared to the local disc appears near the 
planet. As mentionned earlier, this yields a positive feedback on 
migration and possibly gives rise to high drift rates if the viscosity 
and thermal diffusivity are small enough so that the entropy excess 
in the librating region is conserved in the course of migration. 
In the case where the entropy deficit is negative, however, the 
librating region tends to be cooler and over dense relative to the 
ambient disc, resulting in a negative feedback on migration. In that 
case, we expect migration in non-isothermal discs to be strongly 
slowed down. 

As mentionned earlier, we expect this mecanism to arise 
whenever the drift timescale across the horseshoe region is shorter 
than the libration time. This condition reads (Papaloizou et al. 
2007): 


a > df 


(16) 


with (Peplinski et al. 2008); 


df 


Sna 


(17) 


Using Eq. [T^ we can now estimate under which condition on the 
disc mass, horseshoe streamlines can be significantly distorted by 
the effect of the radial drift of the planet. The migration rate d of 
the planet is (e.g. Paardekooper 2014): 


d = la — 

trip 


F/Fo ^ 

, i 0 
yJGM^,a 


(18) 


Given that Fo = {qjhpf-'Lpa'^tFp, the previous equation can be re¬ 
written as: 


a = -(F/Fo)^flnp (19) 

n hi 

Substituting the latter expression in Eq. M and using x, ~ 
l-la-fq/lf, we find that a protoplanet may enter a fast migration 
regime provided that: 


Si > 

hp (yF/Fo)’ 


( 20 ) 




Figure 4. Upper panel. Disc aspect ratio as a function of radius for the 
two radiative disc models that are considered and with parameters given in 
Table 1. Lower panel. Toomre stability pai'ameter as a function of radius for 
the two models. The dashed line corresponds to the gravitational stability 
limit. 


or equivalently: 


Q < —(yF/Fo), (21) 

y 

where Q = hlqj is the Toomre parameter, and where yF/Fo is given 
by Eq.|^ For y = 1.4 and using p = ^ + {y — l)cr, Eq.|^reads: 


yF/Fo = -2.5 -H 3.9f - 0.6cr (22) 

Not surprisingly, this predicts that it is easier for the planet to enter 
the fast migration regime i) for high values of the entropy gradient, 
and ii) when the planet migrates inward and the entropy gradient 
is positive. This also shows that in the case where the planet mi¬ 
grates in the opposite direction to that predicted by the corotation 
torque, namely if ^ has a moderate and positive value, yF/Fo is of 
the order of unity and distortion of the horseshoe streamlines oc¬ 
curs only provided that the disc is relatively massive. For the disc 
model with ^ = -1, we note that Eq.|20]predicts fast migration for 
> 0.003 which is in perfect agreement with the results of the 
simulations. In the case with ^ = 2, however, Eq. [^predicts fast 
migration for q^ > 0.005 whereas the right panel of Fig. [T] shows 
fast migration for qj > 0.01. This discrepancy possibly arises be¬ 
cause the expression given by Eq. does not take into account 
possible saturation effects. We therefore expect this estimation to 
slightly underestimate the critical disc mass above which fast mi¬ 
gration should occur. 
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Figure 5. Upper panel: Contour plots showing, for the two radiative disc models, the disc torque as a function of planet mass and planet radial location. This 
has been derived using the analytical formulae of Paardekooper et al. (2011). Lower panel: this shows, for the two radiative discs, the time evolution of the 
semi-major axes of planets of different masses and that are initially located in the region where outward migration is expected to proceed. 


4 RADIATIVE DISC MODELS 
4.1 Results 

We now turn to the case of equilibrium radiative discs for which 
the thermodynamical state is determined by the balance between 
viscous heating, stellar irradiation, and radiative cooling. We plot 
in Fig. |4] the disc aspect ratio (upper panel) and Toomre stability 
parameter (lower panel) as a function of radius for the two mod¬ 
els we consider. We remind the reader that the disc parameters for 
the two models can be found in Table 1. For Model 1, the struc¬ 
ture in the disc aspect ratio is very similar to that of stellar equilib¬ 
rium discs (SED) obtained using 3D hydrodynamical simulations 
(see for example Fig. 1 in Lega et al. 2015). In the inner disc, vis¬ 
cous heating dominates over stellar irradiation and opacity transi¬ 
tions create bumps in HjR si R ~ 0.6 and ~ 1 (Bitsch et al. 
2013, 2014). In the outer parts of the disc, however, stellar irradi¬ 
ation is the main source of heating, leading to a flared disc with 
HjR oc fQj- R > 2. For Model 2, the aspect ratio is slightly 
smaller due to the lower value for the viscosity but the radial pro¬ 
file of HIR looks similar to that for Model 1. For both models, the 
torque as a function of planet mass and orbital distance is presented 
in the upper panel of Fig.|5] Following Bitsch et al. (2013), these 
migration maps have been obtained using the analytical formulae 
of Paardekooper et al. (2011) which include saturation effects for 
the corotation torque. Here, the black lines delimitate the regions 
in the disc where outward migration should occur. As expected, 
Type I migration is directed inward in the flared parts of the disc 


where the entropy gradient is positive, yielding a negative entropy- 
related corotation torque. In the regions where the disc aspect ra¬ 
tio decreases with radius, however, the entropy gradient tends to 
be negative, and this can possibly lead to outward migration if the 
corotation torque is not saturated and the entropy gradient steep 
enough. Because the thermal diffusion timescale tends to increase 
with radius, outward migration will proceed in that case until the 
planet reaches the zero-torque radius where the saturated corota¬ 
tion torque counterbalances the Lindblad torque. As mentionned 
in Sect. 13.11 this occurs when the thermal diffusion timescale is 
approximately equal to the libration timescale in the horseshoe re¬ 
gion. Moreover, since the libration period decreases as the planet 
mass increases, the zero-torque radius set by saturation tends to be 
located closer in for heavier planets, which is clearly seen in the 
upper panel of Fig. |5] For Model 1 (resp. Model 2) , protoplan¬ 
ets with q < 6 X 10“= (resp. q < 1.6 X 10“=) would have their 
zero-torque radius located in the flared part of the disc where the 
entropy gradient is positive. For this mass range, the location of the 
zero-torque radius becomes almost mass-independent (Cossou et 
al. 2013) and rather corresponds to the transition between the vis¬ 
cous heating dominated regime and the stellar heating dominated 
region. This is located st R ~ 2 for Model 1 while it is located at 
I? ~ 1.5 for Model 2. 

We note that these migration maps have been computed 
from analytical formulae that are valid for planets held on a 
fixed circular orbit only. Returning to Fig. |4l however, we see 
that the Toomre stability parameter can reach values in the range 
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2.5 < < 5 in the regions where outward migration is expected 

to occur. From the discussion presented in Sect. 13.21 . it cannot he 
excluded that the criterion given hy Eq.[2T]is fulfilled in this case, 
and that drift rates significantly higher than those predicted hy 
static torques can he reached. More precisely, for both models, the 
entropy gradient in the region of outward migration is estimated 
to be such that ^ ~ 1.8, which gives yF/Fo ~ 4 from Eq.[^ The 
condition given by Eq. ED therefore becomes < 10, which is 
clearly fulfilled for Model 2, and for Model 1 as well in regions 
withf? > 1.3 (see lower panel of Eig.|4j. 

Eor Model 1, we present in the lower left panel of Eig. 

the semi-major axis as a function of time for planets with 
q - 3 X 10“^, 6 X 10“^, 10“'* and initially located at ao = 1.1, 
namely at a location where outward migration is expected to occur 
for this mass range. The orbital evolution of the q = 3 x 10“^ pro¬ 
toplanet is in close agreement with that predicted from the upper 
panel of Eig. involving outward migration until the zero-torque 
radius located at ~ 2.1 is reached. Eor q = 6 x 10“^, 10“'*, 
however, the evolution differs significantly to what expected from 
the migration map. The drift rates for these cases appear to be 
significantly higher than those estimated using static saturated 
or unsaturated corotation torques, leading the planets to migrate 
outward well beyond the zero-torque radius. This is exemplified in 
Eig.|6] where we plot the planet drift rate as a function of radius 
for both cases. Eor the run with q = 10^“*, it is clear that the 
planet enters a very fast migration regime once R ~ 1.8, which 
corresponds to a value for the Toomre parameter of 4 (see the 
lower panel of Eig. |4j. We suspect this value to correspond to the 
limit below which the mechanism presented in Sect. l3.2l is at work. 
To demonstrate that this is indeed the case, we show in Eig. |7] a 
series of surface density plots at different times for this simulation. 
At early times, there is a negative (resp. positive) surface density 
perturbation at the outward (resp. inward) downstream separatrix 
due to the negative entropy gradient, and the horseshoe region is 
clearly visible. At f ~ 350 orbits, namely just before the planet 
starts its fast outward migration, an underdense librating region has 
appeared, suggesting thereby that the condition given by Eg.llhlis 
verified at that time. As mentionned in Sect. 13.21 this can yield a 
strong positive corotation torque which is responsible for the high 
drift rates observed. As can be seen in the third panel of Eig. |7] 
the trapped region progressively shrinks in the course of migration 
and a part of the coorbital mass deficit can be lost, resulting in the 
planet migrating inward again at later times. Interestingly, we note 
that a similar process prevents outward runaway migration (Masset 
& Papaloizou 2003) to be sustained over long timescales. Here, 
thermal diffusion acts also against fast migration since it makes 
the trapped region cool down, leading to a slight increase in the 
surface density. Once the planet migrates inward again and evolves 
in the flared part of the disc, the thermodynamical state of the disc 
is close to the isothermal limit so that radiative effects can not lead 
to the formation of an underdense librating region near the planet, 
resulting in the planet migrating under the action of static torques 
only. 

Eor Model 2, the orbital evolution of planets with q - 
10“^, 1.5 X 10“^, 2 X 10“^ is shown in the lower right panel of Eig. 

Outward migration here occurs for ^ < 1.5 x 10“^ whereas satu¬ 
ration of the corotation torque makes the q = 2x 10“^ protoplanet 
migrate inward. In the latter case, we note that this does not agree 
with the migration map in the upper panel of Eig. which pre¬ 
dicts outward migration for R < 1.2. Inspection of this migration 



Figure 6. Drift rate as a function of the radial location of the planet for the 
simulations with q = 6x 10^^, 10^"* and for the radiative disc corresponding 
to Model 1. The dashed line represents the drift rate that is obtained using 
the formulae of Paardekooper et al. (2011) for the disc torque, whereas the 
dot-dashed line represents the drift rate predicted assuming a fully unsatu¬ 
rated corotation torque. 

map also reveals that compared to the case with q = 1.5 x 10“^, 
the level of saturation is clearly initially smaller for a planet mass 
with q - 10“^. The resulting faster outward migration makes the 
q = 10“^ body rapidly reach regions with relatively low Q values. 
Analysis of the streamlines reveals that the horseshoe region is de¬ 
stroyed at f 700 orbits, at which time the planet has migrated to a 
disc region where Q ~ 3.5. Again, this is fairly consistent with the 
estimation given by Eq. [21] Erom that time, the planet undergoes 
an episode of fast migration due to the effect presented above, but 
becomes ultimately trapped at the zero-torque radius located at the 
transition between the viscous heating dominated and stellar heat¬ 
ing dominated regimes. In agreement with the discussion of Sect. 
13.21 we find that the typical drift rates during fast migration can be 
higher than those predicted assuming a fully unsaturated corotation 
torque, as illustrated by Eig. [8] which shows for this run the drift 
rate of the planet as a function of its radial location in the disc. 

Eor q = 6x 10“^, we mention that continuation of the run in¬ 
dicates that the planet does not appear to experience fast migration, 
in spite of the fact that it can evolve in disc regions with Q ~ 3.5. 
This possibly arises because, as mentionned above, the corotation 
torque is quite saturated in that case, resulting in a low contrast of 
entropy initially between the trapped material and the local disc. 

4.2 Implication for the formation of giant planet cores 

It has been suggested that the zero-torque radius may represent an 
ideal site for the growth of giant planet cores through giant im¬ 
pacts between Earth-mass embryos (Lyra et al. 2010; Hasegawa & 
Pudritz 2011). Pierens et al. (2013) have studied the evolution of 
multiple planets of a few Earth masses which migrate convergently 
toward a convergence zone located at the transition between two 
different opacity regimes. They found that the embryos tend to form 
a stable resonant chain that protect them from close encounters with 
other bodies, preventing thereby the subsequent formation of giant 
planet cores. In the case of a large number of initial embryos or if 
the effects of disc turbulence are included these resonant configura¬ 
tions can be disrupted, but the formation of massive cores remains 
marginal. In a subsequent study, Zhang et al. (2014) focused on 
the case where the zero-torque radius is located at the transition 
between the inner viscously heated part of the disc and the outer 
irradiated region. They found that in massive discs with accretion 
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Figure 7. Maps of the perturbed surface density at f = 120,350,500 orbits for the simulation with ^ = 10 and for the radiative disc coiTesponding to 
Model 1. 


rates M ~ 10 M^lyr the resonant barriers that are formed can be 
broken, enabling the formation of massive cores. 

In massive discs, the mechanism presented in Sect. 13.21 may 
also facilitate the disruption of resonances at convergence zones 
and consequently promote collisions between embryos. To inves¬ 
tigate this issue in more details, we have examined the evolution 
of four planets with q = 10“^ embedded in the disc corresponding 
to Model 2, and located in the vicinity of the zero-torque radius at 
R ~ 1.5. According to the results presented in the previous sec¬ 
tion, the bodies located inside the zero-torque radius are expected 
to rapidly migrate outward. 

Here, we used a setup similar to that presented in Pierens et al. 
(2013), with two embryos initially located on each side of the zero- 
torque radius and separated by 4.5 R,„h, where R„h is the mutual 
Hill radius. The protoplanets are initially set on moderately inclined 
orbits, and a prescription for the rate of inclination damping due to 
the interaction with the disc is employed (see Pierens et al. 2013). 

Fig-Elshows the time evolution of the semi major axis for the 
4 protoplanets. At early times, the two innermost (resp. outermost) 
bodies located inside the zero-torque radius undergo outward (resp. 
inward) migration, as expected. For the two innermost bodies, fast 
outward migration is observed, in agreement with the results of 
the previous section. This causes the embryo initially located at 
rtfl = 1-1 (blue) to rapidly reach the zero-torque radius and to enter 
in a 9 ; 8 resonance with the third planet (black) at r ~ 700 orbits. 
The inward migrating fourth planet (red) then catches up with these 
two bodies and enter in a 9 : 8 resonance with them. This resonant 
chain is however subsequently destabilized due to the fast outward 
migration of the first body (green), resulting in a physical collision 
between the two outermost planets at r ~ 1000 orbits. The 6.6 
planet that is formed at that time by this process subsequently mi¬ 
grates inward until it reaches the zero-torque radius at ~ 1.5. 
Regarding the two innermost embryos, they undergo a series of or¬ 
bital exchanges inside the zero-torque until they collide at f ~ 2800 
orbits, leading again to the formation of a 6.6 planet. From the 
migration map corresponding to model 2, we expect this body to 
slowly migrate inward at later times until it reaches the zero-torque 
radius set by saturation and which is located at R ~ 1.2 for such a 
planet mass. 

This result seems to indicate that in relatively massive discs for 
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Figure 8. Drift rate as a function of the radial location of the planet for 
the simulations with q = 10“^ and for the radiative disc corresponding to 
Model 2. The dashed line represents the drift rate that is obtained using the 
formulae of Paardekooper et al. (2011) for the disc torque, whereas the dot- 
dashed line represents the drift rate predicted assuming a fully unsaturated 
corotation torque. 

which the condition given by Eq. ED is fulfilled, fast convergent 
migration of low-mass embryos may indeed enhance the formation 
of massive cores through collisional growth of low-mass bodies. 
We will address in a future publication the issue of the dependence 
of these results on the disc model and initial mass of the embryos, 
and investigate whether or not giant planet cores may indeed be 
formed by this process. 


5 DISCUSSION AND CONCLUSION 

In this paper, we have presented the results of hydrodynamic sim¬ 
ulations of the orbital evolution of embedded low-mass planets in 
non-isothermal and radiative discs. A main aim of this study is to 
examine the dependence of the planet drift rate on the disc mass 
and entropy gradient inside the disc, with particular emphasis put 
on the possible role of dynamical corotation torques. We consider 
non-barotropic discs with an initial radial entropy gradient, and in 
which both the thermal diffusion coefficient and viscosity are con- 
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Figure 9. Orbital evolution of four embryos with mass ratio g = 10“^ ini¬ 
tially located near the transition between the viscous heating dominated and 
stellar heating dominated regimes, for the radiative disc model correspond¬ 
ing to model 2. 


slant. We also investigate the case of radiative disc models that in¬ 
clude the effect of viscous heating, stellar irradiation and radiative 
cooling. For non-isothermal disc models with an initial non-zero 
entropy gradient and zero vortensity gradient, we find that the drift 
rates of a migrating protoplanet can be, under certain conditions, 
much higher than those expected from classical formulae for the 
disc torque. We observe that these can even be higher than those 
expected assuming a fully unsaturated horseshoe drag. This arises 
provided that i) the Toomre stability parameter is less than a criti¬ 
cal value that depends on the entropy gradient inside the disc and 
which is given by Eq.|^and ii) the planet migrates in the direction 
set by the entropy-related horseshoe drag. If condition i) is satisfied, 
the horseshoe region of the planet does not extend to the full 2n in 
azimuth and rather contracts into a tadpole-like region (Papaloizou 
et al. 2007). This region is bound to the planet and has therefore a 
negative feedback on migration. It is located at the leading side of 
the planet if the latter migrates inward, whereas it is located at trail¬ 
ing side of the planet if migration is inward. However, if condition 
ii) is also verified, the trapped material shows an entropy excess 
and becomes underdense relative to the ambient disc as the disc 
tries to maintain a pressure balance. In that case, the main contri¬ 
bution to the corotation torque comes from the material that flows 
across the orbit as the planet migrates and which yields a positive 
feedback on migration. This can make the planet enter in a fast mi¬ 
gration regime, provided that the viscosity and thermal diffusivity 
are small enough so that the entropy excess in the trapped region 
tends to be conserved in the course of migration. 

As a side result, we find that an outward migrating planet that 
enters such a fast migration mode can pass through the location of 
the zero-torque radius due to saturation effects, in agreement with 
the results of Paardekooper (2014) in isothermal discs. 

Our simulations for radiative disc models essentially confirm 
these findings. For a stellar equilibrium disc with constant viscosity 
V = 10“^ and disc mass equivalent to 5 times the Minimum Mass 
Solar Nebula (MMSN) at 5 AU, planets with mass ratio in the range 
6 X 10“^ < q < 10“"^ are observed to enter a fast migration regime 
in the region 5 < 7? < 10 AU where outward migration should pro¬ 
ceed according to the corresponding migration map (see Fig.j^. In 
agreement with the analytical estimation given by Eq. ED this oc¬ 
curs once the Toomre stability parameter becomes smaller than the 


critical value ~ 3 - 4. From that time, strong corotation torques 
make the planet not only pass through the location of the zero- 
torque radius set by saturation, but also migrate well beyond the 
boundary corresponding to the transition between the viscous heat¬ 
ing dominated regime and the stellar heating dominated regime. For 
example, we observe that fast outward migration of Neptune-mass 
planets can proceed until the planet reaches R ~ 20 AU , whereas 
migration contours predict inward migration outside 7 AU. For a 
model with v = 10“® and disc mass equivalent to 2.5 times the 
MMSN at 5 AU, we find that the mechanism presented above can 
strongly affect the outward migration of protoplanets with q ~ 10“^ 
(equivalent to 3.3 M®). In that case, it appears that the drift rate is 
again higher than the one expected assuming a fully unsaturated 
corotation torque, but the planet becomes ultimately trapped at the 
transition between the viscous heating dominated and stellar heat¬ 
ing dominated regimes. For q - 10“^, we have examined the impact 
of this fast migration regime on the orbital evolution of multiple 
planets that convergently migrate at this transition. We find that the 
high drift rates that are achieved in the course of outward migration 
may help in disrupting the resonant chains that are established be¬ 
tween embryos. The consequence is that the process of giant core 
growth by collisions of low-mass bodies may be easier if these em¬ 
bryos undergo fast migration, but this needs to be investigated in 
more details. 

A major simplification of this work resides in the fact that we 
have not considered the effect disc self-gravity. However, we have 
seen that the typical value for the Toomre stability parameter below 
which fast migration is triggered is typically ~ 3 - 4 in radiative 
discs. For such a value, the shift in the positions of the Lindblad res¬ 
onances caused by self-gravity (Pierens & Hure 2005) is expected 
to lead to a Lindblad torque which is ~ 15% stronger in compari¬ 
son with the case where self-gravity is discarded (Baruteau & Mas- 
set 2008). In the case where the planet migrates outward (resp. in¬ 
ward), we therefore expect Eq. |2T]to slightly underestimate (resp. 
overestimate) the value for Q below which the horseshoe stream¬ 
lines become significantly distorted by the radial drift of the planet. 
Clearly, additional self-gravitating simulations are required to defi¬ 
nitely assess how these results are impacted by including the effect 
of the disc self-gravity. This issue will be examined in a future pub¬ 
lication. 
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